#include <iostream>
#include <stdexcept>
#include <vector>
#include <ginac/ginac.h>
#include "sector_decomposition/sector_decomposition.h"

int main()
{

  try{

    using namespace sector_decomposition;
    using namespace GiNaC;

    std::cout << "testsuite sector decomposition" << std::endl;
    std::cout << std::endl;

    // --------------------------------------------------------------
    int verbose_level = 1;

    CHOICE_STRATEGY = STRATEGY_C;

    // --------------------------------------------------------------
    monte_carlo_parameters mc_parameters = monte_carlo_parameters( 5, 15, 10000, 100000 );

    // --------------------------------------------------------------

    symbol eps("eps");
    int order = -2;

    int n     = 3;
    int loops = 1;

    symbol minus_p12("minus_p12"), minus_p22("minus_p22"), minus_p32("minus_p32"), m12("m12"), m22("m22"), m32("m32");

    std::vector<ex> parameters;
    symbol x1("x1"), x2("x2"), x3("x3");
    parameters.push_back(x1);
    parameters.push_back(x2);
    parameters.push_back(x3);

    integration_data global_data = integration_data(parameters, eps, order);

    std::vector<exponent> nu_minus_1(n,exponent(0,0));

    ex U = x1+x2+x3;
    ex F = x1*x2* minus_p32 + x2*x3* minus_p12 + x3*x1* minus_p22 
      + x1*m12 + x2*m22 + x3*m32;

    F = F.subs( minus_p12 == 1 );
    F = F.subs( minus_p22 == 0 );
    F = F.subs( minus_p32 == 0 );
    F = F.subs( m12 == 0 );
    F = F.subs( m22 == 0 );
    F = F.subs( m32 == 0 );

    std::vector<ex> poly_list;
    poly_list.push_back(U);
    poly_list.push_back(F);

    std::vector<exponent> c(poly_list.size());
    c[0] = exponent( n-(loops+1)*2, loops+1 );
    c[1] = exponent( -n+2*loops, -loops );
    for (int k=0; k<n; k++) 
      {
	c[0].sum_up(nu_minus_1[k]);
	c[1].subtract_off(nu_minus_1[k]);
      }

    integrand my_integrand = integrand(nu_minus_1, poly_list, c);

    // --------------------------------------------------------------

    monte_carlo_result res = do_sector_decomposition(global_data, my_integrand, mc_parameters, verbose_level);

    // --------------------------------------------------------------

  } catch (std::exception &e)
    {
      std::cout << "Exception : " << e.what() << std::endl;
    }

  return 0;
}

